package scu.maqiang.mesh;

import scu.maqiang.numeric.CG;
import scu.maqiang.numeric.Direct;
import scu.maqiang.numeric.MVO;
import scu.maqiang.numeric.ParamCheck;

import java.io.File;
import java.io.FileNotFoundException;
import java.io.FileWriter;
import java.io.IOException;
import java.util.ArrayList;
import java.util.Scanner;

import static java.lang.Math.*;

public class MeshGenerator2D {
    public Mesh2T3 readMesh(String fileName, MeshFormat mf) {
        Mesh2T3 mesh = new Mesh2T3();
        Scanner sc = null;
        String line = "";
        int nv, nt, nb;
        try {
            sc = new Scanner(new File(fileName));
        } catch (FileNotFoundException e) {
            e.printStackTrace();
        }

        mesh.nDim = 2;
        mesh.nPerElement = 3;
        mesh.nPerBoundary = 2;
        mesh.nBoundaryPerElement = 3;

        switch (mf) {
            case msh:
                nv = sc.nextInt();
                nt = sc.nextInt();
                nb = sc.nextInt();
//			mesh.nv = nv;
//			mesh.nt = nt;
//			mesh.nb = nb;
//			mesh.nodes = new double[nv][];
//			mesh.nodeLabel = new int[nv];
//			mesh.elements = new int[nt][];
//			mesh.elementLabel = new int[nt];
//			mesh.boundarys = new int[nb][];
//			mesh.boundaryLabel = new int[nb];
                mesh.initMesh0(nv, nt, nb);
                for (int i = 0; i < nv; i++) {
                    double x = sc.nextDouble();
                    double y = sc.nextDouble();
                    mesh.nodes[i] = new double[] {x, y};
                    mesh.nodeLabel[i] = sc.nextInt();
                }

                for (int i = 0; i < nt; i++) {
                    int ti = sc.nextInt() - 1;
                    int tj = sc.nextInt() - 1;
                    int tk = sc.nextInt() - 1;
                    mesh.elements[i] = new int[] {ti, tj, tk};
                    mesh.elementLabel[i]  = sc.nextInt();
                }

                for (int i = 0; i < nb; i++) {
                    int bi = sc.nextInt() - 1;
                    int bj = sc.nextInt() - 1;
                    mesh.boundaries[i] = new int[] {bi, bj};
                    mesh.boundaryLabel[i] = sc.nextInt();
                }

                sc.close();
                break;
            case mmsh:
                nv = sc.nextInt();
                nt = sc.nextInt();
                mesh.nv = nv;
                mesh.nt = nt;
                mesh.nodes = new double[nv][];
                mesh.nodeLabel = new int[nv];
                mesh.elements = new int[nt][];
                mesh.elementLabel = new int[nt];
                for (int i = 0; i < nv; i++) {
                    double x = sc.nextDouble();
                    double y = sc.nextDouble();
                    mesh.nodes[i] = new double[] {x, y};
                    mesh.nodeLabel[i] = sc.nextInt();
                }

                for (int i = 0; i < nt; i++) {
                    int ti = sc.nextInt() - 1;
                    int tj = sc.nextInt() - 1;
                    int tk = sc.nextInt() - 1;
                    mesh.elements[i] = new int[] {ti, tj, tk};
                    mesh.elementLabel[i]  = sc.nextInt();
                }
                mesh.extractBoundarysFromElements();
                sc.close();
                break;
            case Inp:
                do {
                    line = sc.nextLine();
                }while(!line.startsWith("*Node"));
                line = sc.nextLine();
                ArrayList<double[] > nodes = new ArrayList<>();
                while (!line.startsWith("*Element")){
                    String[] nodesCoord = line.split(",");
                    double x = Double.parseDouble(nodesCoord[1].trim());
                    double y = Double.parseDouble(nodesCoord[2].trim());
                    nodes.add(new double[]{x, y});
                    line = sc.nextLine();
                }
                nv = nodes.size();
                mesh.nv = nv;
                mesh.nodes = new double[nv][];
                mesh.nodeLabel = new int[nv];
                int idx = 0;
                for(double[] xy: nodes) {
                    mesh.nodes[idx++] = xy.clone();
                }
                line = sc.nextLine();
                ArrayList<int[]> elements = new ArrayList<>();
                while(!line.startsWith("*Nset")) {
                    String[] connectivity = line.split(",");
                    int t1 = Integer.parseInt(connectivity[1].trim()) - 1;
                    int t2 = Integer.parseInt(connectivity[2].trim()) - 1;
                    int t3 = Integer.parseInt(connectivity[3].trim()) - 1;
                    elements.add(new int[]{t1, t2, t3});
                    line = sc.nextLine();
                }
                nt = elements.size();
                mesh.elements = new int[nt][];
                idx = 0;
                for(int[] ele: elements) {
                    mesh.elements[idx++] = ele.clone();
                }
                mesh.elementLabel = new int[nt];
                mesh.extractBoundarysFromElements();
                //System.out.println(mesh.nv + "\t" + mesh.nt + "\t" + mesh.nb);
            default:
                break;
        }
        return mesh;
    }

    public Mesh2T3 readMesh(String xyFile, String eleFile) {
        Mesh2T3 mesh = new Mesh2T3();
        Scanner scXY = null;
        Scanner scEle = null;
        FileWriter fwOther = null;
        String line = "";
        try {
            scXY = new Scanner(new File(xyFile));
            scEle = new Scanner(new File(eleFile));
        } catch (IOException e) {
            e.printStackTrace();
        }

        ArrayList<double[]> nodes = new ArrayList<>();
        while(scXY.hasNext()) {
            line = scXY.nextLine();
            String[] xy = line.split(",");
            double t1 = Double.parseDouble(xy[0]);
            double t2 = Double.parseDouble(xy[1]);
            nodes.add(new double[] {t1, t2});
        }
        ArrayList<int[]> eles = new ArrayList<>();
        while(scEle.hasNext()) {
            line = scEle.nextLine();
            String[] ele = line.split(",");
            int t1 = Integer.parseInt(ele[0]);
            int t2 = Integer.parseInt(ele[1]);
            int t3 = Integer.parseInt(ele[2]);
            eles.add(new int[] {t1, t2, t3});
        }
        int nv = nodes.size();
        int nt = nodes.size();
        int nb = 0;
        mesh.initMesh0(nv, nt, nb);
        for(int i = 0; i < nv; i++) {
            mesh.nodes[i] = nodes.get(i);
        }
        for(int i = 0; i < nt; i++) {
            mesh.elements[i] = eles.get(i);
        }
        mesh.extractBoundarysFromElements();

        return mesh;
    }

    public Mesh2T3 square2D(int xNum, int yNum, T3Type type) {
        Mesh2T3 mesh = new Mesh2T3();
        int nv = (xNum + 1) * (yNum + 1);
        int nt = 2 * xNum * yNum;
        int nb = 2 * (xNum + yNum);
        mesh.nv = nv;
        mesh.nt = nt;
        mesh.nb = nb;

        mesh.initMesh0(nv, nt, nb);
        int count = 0;
        for (int i = 0; i < yNum + 1; i++) {
            for (int j = 0; j < xNum + 1; j++) {
                mesh.nodes[count++] = new double[] {(j + 0.0) / xNum,(i + 0.0) / yNum};
            }
        }

        createSquareElements(mesh, xNum, yNum, type);
        createSquareBoundarys(mesh, xNum, yNum, Direct.X);
        return mesh;
    }

    public Mesh2T3 square2D(int xNum, int yNum, T3Type type, Direct direct) {
        if (direct == Direct.X) {
            return square2D(xNum, yNum, type);
        } else if (direct == Direct.Y) {
            Mesh2T3 mesh = new Mesh2T3();
            int nv = (xNum + 1) * (yNum + 1);
            int nt = 2 * xNum * yNum;
            int nb = 2 * (xNum + yNum);
            mesh.nv = nv;
            mesh.nt = nt;
            mesh.nb = nb;
            int count = 0;
            mesh.initMesh0(nv, nt, nb);
            for (int j = 0; j < xNum + 1; j++) {
                for (int i = 0; i < yNum + 1; i++) {
                    mesh.nodes[count++] = new double[] {(j + 0.0) / xNum, (i + 0.0) / yNum};
                }
            }
            createSquareElements(mesh, yNum, xNum, type, direct);
            createSquareBoundarys(mesh, yNum, xNum, direct);
            return mesh;
        } else {
            throw new IllegalArgumentException("The parameter Direct must only be Direct.X or Direct.Y!");
        }
    }

    public Mesh2T3 square(int xNum, int yNum) {
        int t1, t2, t3, t4, t5, nv, nt, nb;
        nv = 2 * xNum * yNum + xNum + yNum + 1;
        nt = 4 * xNum * yNum;
        nb = 2 * (xNum + yNum);
        Mesh2T3 mesh = new Mesh2T3();
        mesh.initMesh0(nv, nt, nb);

        int count = 0;
        for(int j = 0; j < yNum; j++) {
            for(int i = 0; i < xNum + 1; i++) {
                mesh.nodes[count++] = new double[] {(i + 0.0) / xNum, (j + 0.0) / yNum};
            }
            for(int i = 0; i < xNum; i++) {
                mesh.nodes[count++] = new double[] {(i + 0.5) / xNum, (j + 0.5) / yNum};
            }
        }

        for(int i = 0; i < xNum + 1; i++) {
            mesh.nodes[count++] = new double[]{(i + 0.0) / xNum, 1.0};
        }

        count = 0;
        for(int j = 0; j < yNum; j++) {
            for(int i = 0; i < xNum; i++) {
                t1 = i + j * (2 * xNum + 1);
                t2 = t1 + 1;
                t3 = t2 + 2 * xNum + 1;
                t4 = t3 - 1;
                t5 = t1 + xNum + 1;
                mesh.elements[count++] = new int[] {t1, t2, t5};
                mesh.elements[count++] = new int[] {t2, t3, t5};
                mesh.elements[count++] = new int[] {t3, t4, t5};
                mesh.elements[count++] = new int[] {t4, t1, t5};
            }
        }

        count  = 0;
        for(int i = 0; i < xNum; i++) {
            mesh.boundaries[count] = new int[] {i, i + 1};
            mesh.boundaryLabel[count++] = 1;
        }

        int local = xNum;
        for(int i = 0; i < yNum; i++) {
            mesh.boundaries[count] = new int[] {local, local + 2 * xNum + 1};
            mesh.boundaryLabel[count++] = 2;
            local += 2 * xNum + 1;
        }

        local = 2 * xNum * yNum + xNum + yNum;
        for(int i = 0; i < xNum; i++) {
            mesh.boundaries[count] = new int[] {local, local - 1};
            mesh.boundaryLabel[count++] = 3;
            local--;
        }



        local = yNum * (2 * xNum + 1);
        for(int i = 0; i < yNum; i++) {
            mesh.boundaries[count] = new int[] {local, local - 2 * xNum - 1};
            mesh.boundaryLabel[count++] = 4;
            local -= 2 * xNum + 1;
        }
        return mesh;
    }

    public Mesh2T3 quadrilateral(double[] QuadX, double[] QuadY, int xNum, int yNum) {
        Mesh2T3 mesh = new Mesh2T3();

        //锟斤拷锟斤拷锟斤拷锟斤拷锟斤拷

        int n, i, j, xNodes, yNodes, ti, tj, tk, tl;
        double xi, et, x, y;
        double[] h = new double[4];
        double[] dXi = new double[4];
        double[] dEt = new double[4];
        double[] hh = new double[8];
        double[] ddXi = new double[8];
        double[] ddEt = new double[8];
        double[] tji = new double[2];
        double[] tli = new double[2];
        n = QuadX.length;
        xNodes = xNum + 1;
        yNodes = yNum + 1;
        int nv = xNodes * yNodes;
        int nt = 2 * xNum * yNum;
        int nb = 2 * (xNum + yNum);
        mesh.initMesh0(nv, nt, nb);
        int count = 0;
        if (n == 4) {
            for (j = 0; j < yNodes; j++) {
                for (i = 0; i < xNodes; i++) {
                    xi = 2 * (i + 0.0) / xNum - 1.0;
                    et = 2 * (j + 0.0) / yNum - 1.0;
                    CG.Q4ShapeFunction(xi, et, h, dXi, dEt);
                    x = MVO.dot_product(h, QuadX);
                    y = MVO.dot_product(h, QuadY);
                    mesh.nodes[count++] = new double[] {x, y};
                }
            }
        } else if (n == 8) {
            for (j = 0; j < yNodes; j++) {
                for (i = 0; i < xNodes; i++) {
                    xi = 2 * (i + 0.0) / xNum - 1.0;
                    et = 2 * (j + 0.0) / yNum - 1.0;
                    CG.Q8ShapeFunction(xi, et, hh, ddXi, ddEt);
                    x = MVO.dot_product(hh, QuadX);
                    y = MVO.dot_product(hh, QuadY);
                    mesh.nodes[count++] = new double[] {x, y};
                }
            }
        } else
            throw new IllegalArgumentException("The length of the array must be 4 or 8!");

        count = 0;
        for (j = 0; j < yNum; j++) {
            for (i = 0; i < xNum; i++) {
                ti = i + j * xNodes;
                tj = ti + 1;
                tk = tj + xNodes;
                tl = tk - 1;
                for (int kk = 0; kk < 2; kk++) {
                    tji[kk] = mesh.nodes[tj][kk] - mesh.nodes[ti][kk];
                    tli[kk] = mesh.nodes[tl][kk] - mesh.nodes[ti][kk];
                }
                if (MVO.dot_product(tji, tli) > 0.0) {
                    mesh.elements[count++] = new int[] {ti, tj, tl};
                    mesh.elements[count++] = new int[] {tj, tk, tl};
                } else {
                    mesh.elements[count++] = new int[] {tj, tk, ti};
                    mesh.elements[count++] = new int[] {tk, tl, ti};
                }
            }
        }

        createSquareBoundarys(mesh, xNum, yNum);
        return mesh;
    }

    public Mesh2T3 ring(double ra, double rb, int nR, int nL) {
        Mesh2T3 mesh = new Mesh2T3();
        double[] theta = MVO.linspace(0, 2 * PI, nR + 1);
        int nv = (nL + 1) * nR;
        int nt = 2 * nR * nL;
        int nb = 2 * nR;
        mesh.initMesh0(nv, nt, nb);
        double dr = (rb - ra) / nL;
        double x, y;
        int count = 0;
        for (int i = 0; i < nL + 1; i++) {
            for (int j = 0; j < nR; j++) {
                x = (ra + i * dr) * cos(theta[j]);
                y = (ra + i * dr) * sin(theta[j]);
                // System.out.println(MVector.toString(x, y));
                mesh.nodes[count++] = new double[] {x, y};
            }
        }

        count = 0;
        int L1, L2, L3, L4;
        for (int i = 0; i < nL; i++) {
            for (int j = 0; j < nR - 1; j++) {
                L1 = j + i * nR;
                L2 = L1 + nR;
                L3 = L2 + 1;
                L4 = L1 + 1;
                mesh.elements[count++] = new int[] {L1, L2, L3};
                mesh.elements[count++] = new int[] {L1, L3, L4};
            }
            L1 = i * nR + nR - 1;
            L2 = L1 + nR;
            L3 = (i + 1) * nR;
            L4 = i * nR;
            mesh.elements[count++] = new int[] {L1, L2, L3};
            mesh.elements[count++] = new int[] {L1, L3, L4};
        }

        count = 0;
        for (int i = 0; i < nR - 1; i++) {
            mesh.boundaries[count] = new int[] {i + 1, i};
            mesh.boundaryLabel[count++] = 1;
        }
        mesh.boundaries[count] = new int[] {0, nR - 1};
        mesh.boundaryLabel[count++] = 1;

        int local = nL * nR;
        for (int i = 0; i < nR - 1; i++) {
            mesh.boundaries[count] = new int[] {local, local + 1};
            mesh.boundaryLabel[count++] = 2;
            local++;
        }

        mesh.boundaries[count] = new int[] {nR * (nL + 1) - 1, nL * nR};
        mesh.boundaryLabel[count] = 2;
        return mesh;
    }

    public Mesh2T3 triangle(double[] x, double[] y, int num) {
        double avX = MVO.average(x);
        double avY = MVO.average(y);
        double[] quadX1 = {x[0], 0.5 * (x[0] + x[1]), avX, 0.5 * (x[0] + x[2])};
        double[] quadY1 = {y[0], 0.5 * (y[0] + y[1]), avY, 0.5 * (y[0] + y[2])};
        Mesh2T3 mesh1 = quadrilateral(quadX1, quadY1, num / 2, num / 2);

        double[] quadX2 = {x[1], 0.5 * (x[1] + x[2]), avX, 0.5 * (x[0] + x[1])};
        double[] quadY2 = {y[1], 0.5 * (y[1] + y[2]), avY, 0.5 * (y[0] + y[1])};
        mesh1.add(quadrilateral(quadX2, quadY2, num/2, num/2));

        double[] quadX3 = {x[2], 0.5 * (x[2] + x[0]), avX, 0.5 * (x[1] + x[2])};
        double[] quadY3 = {y[2], 0.5 * (y[2] + y[0]), avY, 0.5 * (y[1] + y[2])};
        mesh1.add(quadrilateral(quadX3, quadY3, num/2, num/2));
        return mesh1;
    }
//	public static Mesh2DT3 combine(Mesh2DT3... mesh) {
//		int n = mesh.length;
//		Mesh2DT3 finalMesh;
//		finalMesh = addMesh(mesh[0], mesh[1]);
//		for (int i = 2; i < n; i++) {
//			finalMesh = addMesh(finalMesh, mesh[i]);
//		}
//		return finalMesh;
//	}

    /**
     * 锟斤拷锟斤拷锟斤拷锟斤拷喜锟斤拷锟斤拷锟斤拷锟斤拷锟� 要锟斤拷锟斤拷锟斤拷锟斤拷锟节边斤拷锟斤拷锟斤拷锟洁互匹锟斤拷
     *
     * @param meshA
     *            锟斤拷锟斤拷A
     * @param meshB
     *            锟斤拷锟斤拷B
     * @return 锟斤拷锟斤拷锟斤拷之锟斤拷
     */
//	public static Mesh2DT3 addMesh(Mesh2DT3 meshA, Mesh2DT3 meshB) {
//		Mesh2DT3 meshC = new Mesh2DT3();
//
//		double tol = 1.0e-3, xA, yA, xB, yB;
//		int i, j, nvA, nvB, nConnected, n, nBd;
//		int[] idx, idx2;
//		int[] temp = new int[2];
//		int[] newNodesB = new int[meshB.nv];
//		int[] newBdB = new int[meshB.nb];
//		int[] ABd = new int[meshA.nb];
//		Arrays.fill(newNodesB, -1);
//		Arrays.fill(newBdB, -1);
//		Arrays.fill(ABd, -1);
//		nvA = meshA.nv;
//		nvB = meshB.nv;
//		for (i = 0; i < nvA; i++) {
//			double[] node = meshA.nodes[i];
//			xA = node[0];
//			yA = node[1];
//			for (j = 0; j < nvB; j++) {
//				node = meshB.nodes[j];
//				xB = node[0];
//				yB = node[1];
//				if (abs(xA - xB) < tol && abs(yA - yB) < tol) {
//					newNodesB[j] = i;
//					break;
//				}
//			}
//		}
//
//		// System.out.println(newNodesB.length);
//		for (i = 0; i < meshB.nb; i++) {
//			idx = meshB.boundarys[i];
//			// System.out.println(Arrays.toString(idx));
//			temp[0] = newNodesB[idx[0]];
//			temp[1] = newNodesB[idx[1]];
//			if (temp[0] > -1 && temp[1] > -1) {
//				newBdB[i] = 1;
//				for (j = 1; j < meshA.nb; j++) {
//					idx2 = meshA.boundarys[j];
//					if ((idx2[0] == temp[0] && idx2[1] == temp[1]) || (idx2[1] == temp[0] && idx2[1] == temp[0])) {
//						ABd[j] = 1;
//						break;
//					}
//				}
//			}
//		}
//		nConnected = 0;
//		n = nvA;
//		for (i = 0; i < nvB; i++) {
//			// 锟斤拷锟紹锟叫的斤拷锟矫伙拷锟斤拷锟紸锟叫的斤拷锟斤拷锟皆ｏ拷锟斤拷锟斤拷锟斤拷锟斤拷锟斤拷B锟叫的斤拷锟斤拷锟斤拷碌锟斤拷锟斤拷锟斤拷械慕锟斤拷锟斤拷
//			if (newNodesB[i] == -1) {
//				newNodesB[i] = n;
//				n++;
//			} else {
//				// 锟斤拷锟斤拷丫锟斤拷锟皆ｏ拷锟斤拷锟斤拷锟斤拷锟斤拷锟皆的斤拷锟侥革拷锟斤拷
//				nConnected++;
//			}
//
//		}
//
//		nBd = 0;
//		n = meshA.nb;
//		for (i = 0; i < meshB.nb; i++) {
//			if (newBdB[i] == -1) {
//				newBdB[i] = n;
//				n++;
//			} else {
//				nBd++;
//			}
//		}
//
//		meshC.nv = nvA + nvB - nConnected;
//		meshC.nt = meshA.nt + meshB.nt;
//		meshC.nb = meshA.nb + meshB.nb - nBd;
//
//		int nvC = meshC.nv;
//		int ntC = meshC.nt;
//		int nbC = meshC.nb;
//		meshC.nodes = new double[nvC][];
//		meshC.nodeLabel = new int[nvC];
//		meshC.elements = new int[ntC][];
//		meshC.elementLabel = new int[ntC];
//		meshC.boundarys = new int[nbC][];
//		meshC.boundaryLabel = new int[nbC];
//		int count = 0;
//		for (i = 0; i < nvA; i++) {
//			meshC.nodes[count] = meshA.nodes[i];
//			meshC.nodeLabel[count] = meshA.nodeLabel[i];
//			count++;
//		}
//
//		for (i = 0; i < nvB; i++) {
//			if (newNodesB[i] > nvA - 1) {
//				meshC.nodes[count] = meshB.nodes[i];
//				meshC.nodeLabel[count] = meshB.nodeLabel[i];
//				count++;
//			}
//		}
//
//		count = 0;
//		for (i = 0; i < meshA.nt; i++) {
//			meshC.elements[count] = meshA.elements[i];
//			meshC.elementLabel[count] = meshA.elementLabel[i];
//			count++;
//		}
//
//		for (i = 0; i < meshB.nt; i++) {
//			int[] ele = meshB.elements[i];
//			int[] newEle = new int[3];
//			for (int kk = 0; kk < 3; kk++) {
//				newEle[kk] = newNodesB[ele[kk]];
//			}
//			meshC.elements[count] = newEle;
//			meshC.elementLabel[count] = meshB.elementLabel[i];
//			count++;
//		}
//
//		for (i = 0; i < meshA.nb; i++) {
//			int[] li = meshA.boundarys[i];
//			meshC.boundarys[count] = li;
//			if (ABd[i] == 1) {
//				meshC.boundaryLabel[count] = 0;
//			}
//			count++;
//		}
//
//		for (i = 0; i < meshB.nb; i++) {
//			if (newBdB[i] != 1) {
//				int[] li = meshB.boundarys[i];
//				int[] newEle = new int[2];
//				for (int kk = 0; kk < 2; kk++) {
//					newEle[kk] = newNodesB[li[kk]];
//				}
//				meshC.boundarys[count] = newEle;
//				meshC.boundaryLabel[count] = meshB.boundaryLabel[i];
//				count++;
//			}
//		}
//		return meshC;
//	}

    /**
     * 锟皆撅拷锟叫撅拷锟斤拷锟斤拷锟剿结构锟斤拷锟斤拷锟斤拷悖拷锟斤拷锟斤拷涞ピ拷锟斤拷锟斤拷锟斤拷锟斤拷呓锟斤拷锟斤拷锟斤拷锟�
     *
     * @param mesh
     *            要锟斤拷锟缴碉拷锟斤拷锟斤拷
     * @param xNum
     *            x锟斤拷锟斤拷锟斤拷锟斤拷锟斤拷锟�
     * @param yNum
     *            y锟斤拷锟斤拷锟斤拷锟斤拷锟斤拷锟�
     * @param type
     *            锟斤拷锟斤拷锟斤拷锟斤拷锟斤拷锟轿凤拷锟斤拷锟斤拷锟斤拷
     */
    private void createSquareElements(Mesh2T3 mesh, int xNum, int yNum, T3Type type) {
        int t1, t2, t3, t4;
        int count = 0;
        for (int i = 0; i < yNum; i++) {
            for (int j = 0; j < xNum; j++) {
                t1 = i * (xNum + 1) + j;
                t2 = t1 + 1;
                t3 = t1 + xNum + 1;
                t4 = t3 + 1;
                switch (type) {
                    case Left:
                        mesh.elements[count++] = new int[] {t1, t2, t4};
                        mesh.elements[count++] = new int[] {t1, t4, t3};
                        break;
                    case Right:
                        mesh.elements[count++] = new int[] {t1, t2, t3};
                        mesh.elements[count++] = new int[] {t2, t4, t3};
                        break;
                    case Cross1:
                        if ((i + j) % 2 == 0) {
                            mesh.elements[count++] = new int[] {t1, t2, t4};
                            mesh.elements[count++] = new int[] {t1, t4, t3};
                        } else {
                            mesh.elements[count++] = new int[] {t1, t2, t3};
                            mesh.elements[count++] = new int[] {t2, t4, t3};
                        }
                        break;
                    case Cross2:
                        if ((i + j) % 2 == 1) {
                            mesh.elements[count++] = new int[] {t1, t2, t4};
                            mesh.elements[count++] = new int[] {t1, t4, t3};
                        } else {
                            mesh.elements[count++] = new int[] {t1, t2, t3};
                            mesh.elements[count++] = new int[] {t2, t4, t3};
                        }
                        break;
                    default:
                        throw new IllegalArgumentException("Unexpected T3Type!");
                }
            }
        }
    }

    private void createSquareElements(Mesh2T3 mesh, int xNum, int yNum, T3Type type, Direct direct) {
        if (direct == Direct.X) {
            createSquareElements(mesh, xNum, yNum, type);
        } else if (direct == Direct.Y) {
            int count = 0;
            for (int i = 0; i < yNum; i++) {
                for (int j = 0; j < xNum; j++) {
                    int ti1, tj1, tk1, ti2, tj2, tk2;
                    switch (type) {
                        case Left:
                            ti1 = i * (xNum + 1) + j;
                            tj1 = ti1 + 1;
                            tk1 = tj1 + xNum + 1;
                            ti2 = i * (xNum + 1) + j;
                            tj2 = ti2 + xNum + 1 + 1;
                            tk2 = tj2 - 1;
                            break;
                        case Right:
                            ti1 = i * (xNum + 1) + j;
                            tj1 = ti1 + 1;
                            tk1 = ti1 + xNum + 1;
                            ti2 = i * (xNum + 1) + j + 1;
                            tj2 = ti2 + xNum + 1;
                            tk2 = tj2 - 1;
                            break;
                        case Cross1:
                            if ((i + j) % 2 == 0) {
                                ti1 = i * (xNum + 1) + j;
                                tj1 = ti1 + 1;
                                tk1 = ti1 + xNum + 1;
                                ti2 = i * (xNum + 1) + j + 1;
                                tj2 = ti2 + xNum + 1;
                                tk2 = tj2 - 1;
                            } else {
                                ti1 = i * (xNum + 1) + j;
                                tj1 = ti1 + 1;
                                tk1 = tj1 + xNum + 1;
                                ti2 = i * (xNum + 1) + j;
                                tj2 = ti2 + xNum + 1 + 1;
                                tk2 = tj2 - 1;
                            }
                            break;
                        case Cross2:
                            if ((i + j) % 2 == 1) {
                                ti1 = i * (xNum + 1) + j;
                                tj1 = ti1 + 1;
                                tk1 = ti1 + xNum + 1;
                                ti2 = i * (xNum + 1) + j + 1;
                                tj2 = ti2 + xNum + 1;
                                tk2 = tj2 - 1;
                            } else {
                                ti1 = i * (xNum + 1) + j;
                                tj1 = ti1 + 1;
                                tk1 = tj1 + xNum + 1;
                                ti2 = i * (xNum + 1) + j;
                                tj2 = ti2 + xNum + 1 + 1;
                                tk2 = tj2 - 1;
                            }
                            break;
                        default:
                            throw new IllegalArgumentException("Unexpected T3Type!");
                    }
                    mesh.elements[count++] = new int[] {ti1, tk1, tj1};
                    mesh.elements[count++] = new int[] {ti2, tk2, tj2};
                }
            }
        } else {
            throw new IllegalArgumentException("The parameter Direct must only be Direct.X or Direct.Y!");
        }
    }

    private void createSquareBoundarys(Mesh2T3 mesh, int xNum, int yNum) {
        int count = 0;
        for (int i = 0; i < xNum; i++) {
            mesh.boundaries[count] = new int[] {i, i + 1};
            mesh.boundaryLabel[count++] = 1;
        }
        int local = xNum;
        for (int i = 0; i < yNum; i++) {
            mesh.boundaries[count] = new int[] {local, local + xNum + 1};
            mesh.boundaryLabel[count++] = 2;
            local += xNum + 1;
        }

        local = xNum + yNum * (xNum + 1);
        for (int i = 0; i < xNum; i++) {
            mesh.boundaries[count] = new int[] {local, local - 1};
            mesh.boundaryLabel[count++] = 3;
            local--;
        }
        local = yNum * (xNum + 1);
        for (int i = 0; i < yNum; i++) {
            mesh.boundaries[count] = new int[] {local, local - (xNum + 1)};
            mesh.boundaryLabel[count++] = 4;
            local -= xNum + 1;
        }
    }

    private void createSquareBoundarys(Mesh2T3 mesh, int xNum, int yNum, Direct direct) {
        if (direct == Direct.X) {
            createSquareBoundarys(mesh, xNum, yNum);
        } else if (direct == Direct.Y) {
            int count = 0;
            for (int i = 0; i < xNum; i++) {
                mesh.boundaries[count] = new int[] {i + 1, i};
                mesh.boundaryLabel[count++] = 4;
            }
            int local = xNum;
            for (int i = 0; i < yNum; i++) {
                mesh.boundaries[count] = new int[] {local + xNum + 1, local};
                mesh.boundaryLabel[count++] = 3;
                local += xNum + 1;
            }

            local = xNum + yNum * (xNum + 1);
            for (int i = 0; i < xNum; i++) {
                mesh.boundaries[count] = new int[] {local - 1, local};
                mesh.boundaryLabel[count++] = 2;
                local--;
            }

            local = yNum * (xNum + 1);
            for (int i = 0; i < yNum; i++) {
                mesh.boundaries[count] = new int[] {local - (xNum + 1), local};
                mesh.boundaryLabel[count++] = 1;
                local -= xNum + 1;
            }
        } else {
            throw new IllegalArgumentException("The parameter Direct must only be Direct.X or Direct.Y!");
        }

    }

    public Mesh2T3 fromRectangleNodes(double[][] nodes, int xNum, int yNum, T3Type type) {
        Mesh2T3 mesh = new Mesh2T3();
        int nv = nodes.length;
        ParamCheck.checkEqual(nv, (xNum + 1) * (yNum + 1));
        int nt = 2 * xNum * yNum;
        int nb = 2 * (xNum + yNum);
        mesh.nv = nv;
        mesh.nt = nt;
        mesh.nb = nb;

        mesh.initMesh0(nv, nt, nb);
        for(int i = 0; i < nv; i++) {
            mesh.nodes[i] = nodes[i].clone();
        }

        createSquareElements(mesh, xNum, yNum, type);
        createSquareBoundarys(mesh, xNum, yNum, xNum <= yNum? Direct.X : Direct.Y);
        return mesh;
    }
}
